第三十讲 R-Cox比例风险模型的假设检验条件

Original 跟我学投必得医学 2022-05-07

收录于合集 #R与生物统计专题 56个

在“R与生物统计专题”中，我们会从介绍R的基本知识展开到生物统计原理及其在R中的实现。以从浅入深，层层递进的形式在投必得医学公众号更新。

在介绍完 Cox比例风险模型的详细理论和R实现以后（第二十八讲 R-Cox比例风险模型（1）、第二十九讲 R-Cox比例风险模型（2）），我们已经知道Cox比例风险模型可以较好的同时矫正多个混杂因素对结果的影响，同时我们遗留下来了一个问题：Cox比例风险模型的假设检验条件是什么？这一讲，我们将对此进行解答。

请注意

如果使用不当，统计模型可能会产生误导性的结论。因此，检查给定的模型是否适合数据类型是非常重要的！

1. Cox比例风险模型的假设条件

Cox比例风险模型的成立是建立在数据满足一定假设条件的基础上的。因此，评估拟合的Cox回归模型是否满足这些假设条件非常重要。

我们需要测试Cox模型成立的三个基本条件：

测试比例风险假设。
检查影响观察结果的值（或异常值）。
检测对数风险值与协变量之间关系的非线性情况。

为了测试这些模型假设，我们推荐使用残差法。

Cox模型的常见残差包括：

Schoenfeld残差：用以检查比例风险假设
Martingale残差：应以评估非线性情况
Deviance残差（Martinguale残差的对称变换）：用以检查异常值

2. 评估 Cox比例风险模型的有效性

2.1 安装和加载所需的R包

我们将使用两个R包：

survival用于计算生存分析
survminer用于可视化生存分析的结果
安装软件包

install.packages(c("survival", "survminer"))

加载包

library("survival")library("survminer")

2.2 计算Cox模型

我们将使用在survival包中的自带lung部数据集和coxph（）函数。

计算Cox模型，请输入以下内容：

library("survival")res.cox <- coxph(Surv(time, status) ~ age + sex + wt.loss, data = lung)res.cox

输出结果

Call:coxph(formula = Surv(time, status) ~ age + sex + wt.loss, data = lung)coef exp(coef) se(coef) z page 0.02009 1.02029 0.00966 2.08 0.0377sex -0.52103 0.59391 0.17435 -2.99 0.0028wt.loss 0.00076 1.00076 0.00619 0.12 0.9024Likelihood ratio test=14.7 on 3 df, p=0.00212n= 214, number of events= 152(14 observations deleted due to missingness)

2.3 测试比例风险假设

我们可以使用一定的统计学方法和图形来描绘标准化（scaled）后的Schoenfeld残差，从而检查比例风险假设。

注意

原则上，Schoenfeld残差与时间无关。如果结果显示该残差与时间有非随机性关系，则说明表明违反了比例风险假设。

函数cox.zph（）[在survival软件包中]提供了一种简便的解决方案，可以对Cox模型拟合中包含的每个协变量进行比例风险假设检验。

对于每个协变量，函数cox.zph（）将标准化的Schoenfeld残差的相应集合与时间相关联，以测试残差和时间之间的独立性。此外，它还会对整个模型进行全局测试。

注意

当比例风险假设中，残差和时间之间的线性关系不显著时，可以说，模型符合比例风险假设。

针对上一步骤得到的比例风险模型，我们进行比例风险假设检验时，请键入以下内容：

test.ph <- cox.zph(res.cox)test.ph

输出结果

rho chisq page -0.0483 0.378 0.538sex 0.1265 2.349 0.125wt.loss 0.0126 0.024 0.877GLOBAL NA 2.846 0.416

从上面的输出中，每个协变量都不具有统计显着性（p > 0.05），全局检验也不具有统计学显着性。因此，我们可以认为该Cox模型符合比例风险假设。

可以使用ggcoxzph（）[在survminer软件包中] 函数进行图形诊断，该函数针对每个协变量生成标准化的Schoenfeld残差相对于时间的相关性图像。

ggcoxzph(test.ph)

在上图中，实线是曲线的平滑样条拟合，虚线代表拟合周围的+/- 2倍标准误差范围。

请注意

如图所示，当线条基本平行于X轴时，我们可以认为符合风险比例假设，因为我们假设的是，系数β1，β2，β3不随时间变化。

从图形检查来看，没有协变量是随时间规律变化的。所有的协变量都满足比例风险的假设。其中，由于性别是一个双数值变量，有男女两个数值，所以图中可以看到有两个条带，分别代表男女。

另一种图形方法来检查比例风险的是绘制log(-log(S(t)))对t或log（t）的相关图。这只能用于分类协变量。其原理类似于上面的方法，这里不做赘述。

那么，如果数据不符合比例风险假设，我们可以通过以下方法解决：

添加协变量*时间交互
分层

分层对于随时间变化的混杂因素很有用，因为这可以有效的避免了这个混杂因素对结果的影响。但是这也使得你没有办法研究这个因素本身对结果的影响。

2.4 检查异常值

为了测试对结果产生影响的异常观察值或离群值，我们可以用下面两个方法可视化：

deviance残差或
dfbeta值

ggcoxdiagnostics（）函数（在survminer软件包中）为检查异常值提供了一种方便的解决方案。简化格式如下：

ggcoxdiagnostics(fit, type = , linear.predictions = TRUE)

fit：类coxph.object的对象
type：要在Y轴上显示的残差类型。允许的值包括c（“ martingale”，“ deviance”，“ score”，“ schoenfeld”，“ dfbeta”，“ dfbetas”，“ scaledsch”，“ partial”）。
linear.predictions：一个逻辑值，指示是显示观测值的线性预测（TRUE）还是仅在X轴上显示观测值的索引（FALSE）。

指定参数类型=“ dfbeta”时，将依次删除每个观察值后绘制回归系数的估计变化；同样，type =“ dfbetas”会产生系数的估计变化除以其标准误。

例如：

ggcoxdiagnostics(res.cox, type = "dfbeta",linear.predictions = FALSE, ggtheme = theme_bw())

（dfbeta的指数图，显示了在Cox回归中，用于死亡时间对年龄，性别和体重损失的影响）

上面的指数图表明，将最大dfbeta值的幅度与回归系数进行比较表明，即使年龄和体重损失的某些dfbeta值相对于其他值而言很大，但这些观察结果与每一个排它结果都没有很大的区别，都是均匀分布在x = 0的参考线两边。

也可以通过可视化deviance残差来检查离群值。deviance残差是martingale残差的正态分布变换形式。当这些残差应大致对称地分布在零附近，标准差应该为1。